{
"cells": [
{
"cell_type": "markdown",
"id": "0",
"metadata": {},
"source": [
"# Lambert Targeting"
]
},
{
"cell_type": "markdown",
"id": "1",
"metadata": {},
"source": [
"This tutorial demonstrates how to use satkit's Lambert solver to design orbit transfers.\n",
"We'll compute the delta-v required for a satellite to move between orbits and visualize the transfer geometry."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "2",
"metadata": {},
"outputs": [],
"source": [
"import matplotlib.pyplot as plt\n",
"import scienceplots # noqa: F401\n",
"plt.style.use([\"science\", \"no-latex\", \"../satkit.mplstyle\"])\n",
"%config InlineBackend.figure_formats = ['svg']\n",
"\n",
"import satkit as sk\n",
"import numpy as np"
]
},
{
"cell_type": "markdown",
"id": "3",
"metadata": {},
"source": [
"## Scenario: LEO to MEO Transfer\n",
"\n",
"A satellite is in a 400 km circular orbit (LEO) and needs to transfer to a 2000 km circular orbit (MEO).\n",
"We'll compute the optimal transfer for a range of transfer angles and find the one that minimizes total delta-v."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "4",
"metadata": {},
"outputs": [],
"source": [
"# Orbit parameters\n",
"alt_1 = 400e3 # departure altitude (m)\n",
"alt_2 = 2000e3 # arrival altitude (m)\n",
"r_earth = sk.consts.earth_radius\n",
"mu = sk.consts.mu_earth\n",
"\n",
"r1_mag = r_earth + alt_1\n",
"r2_mag = r_earth + alt_2\n",
"\n",
"# Circular velocities at each orbit\n",
"v_circ_1 = np.sqrt(mu / r1_mag)\n",
"v_circ_2 = np.sqrt(mu / r2_mag)\n",
"\n",
"print(f\"Departure orbit: {alt_1/1e3:.0f} km altitude, v = {v_circ_1:.1f} m/s\")\n",
"print(f\"Arrival orbit: {alt_2/1e3:.0f} km altitude, v = {v_circ_2:.1f} m/s\")"
]
},
{
"cell_type": "markdown",
"id": "5",
"metadata": {},
"source": [
"## Delta-v vs Transfer Angle\n",
"\n",
"We place the satellite at a fixed departure point and sweep the arrival point around the higher orbit,\n",
"computing the Lambert solution and delta-v for each transfer angle."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "6",
"metadata": {},
"outputs": [],
"source": [
"# Departure position (on x-axis)\n",
"r1 = np.array([r1_mag, 0.0, 0.0])\n",
"\n",
"# Sweep transfer angles from 30 to 170 degrees\n",
"angles_deg = np.linspace(30, 180, 200)\n",
"angles_rad = np.radians(angles_deg)\n",
"\n",
"# Hohmann transfer time (reference)\n",
"a_transfer = (r1_mag + r2_mag) / 2\n",
"t_hohmann = np.pi * np.sqrt(a_transfer**3 / mu)\n",
"\n",
"dv1_arr = []\n",
"dv2_arr = []\n",
"dv_total_arr = []\n",
"\n",
"for theta in angles_rad:\n",
" r2 = np.array([r2_mag * np.cos(theta), r2_mag * np.sin(theta), 0.0])\n",
"\n",
" # Scale TOF with transfer angle: fraction of Hohmann time\n",
" tof = t_hohmann * (theta / np.pi)\n",
"\n",
" solutions = sk.lambert(r1, r2, tof)\n",
" v1_transfer, v2_transfer = solutions[0]\n",
"\n",
" # Departure: satellite in circular orbit, moving in +y direction\n",
" v1_circular = np.array([0.0, v_circ_1, 0.0])\n",
"\n",
" # Arrival: circular velocity tangent to orbit at arrival point\n",
" v2_circular = np.array([\n",
" -v_circ_2 * np.sin(theta),\n",
" v_circ_2 * np.cos(theta),\n",
" 0.0,\n",
" ])\n",
"\n",
" dv1 = np.linalg.norm(v1_transfer - v1_circular)\n",
" dv2 = np.linalg.norm(v2_transfer - v2_circular)\n",
"\n",
" dv1_arr.append(dv1)\n",
" dv2_arr.append(dv2)\n",
" dv_total_arr.append(dv1 + dv2)\n",
"\n",
"dv1_arr = np.array(dv1_arr)\n",
"dv2_arr = np.array(dv2_arr)\n",
"dv_total_arr = np.array(dv_total_arr)\n",
"\n",
"# Find minimum total delta-v\n",
"idx_min = np.argmin(dv_total_arr)\n",
"print(f\"Minimum total delta-v: {dv_total_arr[idx_min]:.1f} m/s at {angles_deg[idx_min]:.1f} degrees\")\n",
"print(f\" Departure burn: {dv1_arr[idx_min]:.1f} m/s\")\n",
"print(f\" Arrival burn: {dv2_arr[idx_min]:.1f} m/s\")"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "7",
"metadata": {},
"outputs": [],
"source": [
"fig, ax = plt.subplots(figsize=(8, 5))\n",
"ax.plot(angles_deg, dv1_arr, label=r\"$\\Delta v_1$ (departure)\")\n",
"ax.plot(angles_deg, dv2_arr, label=r\"$\\Delta v_2$ (arrival)\")\n",
"ax.plot(angles_deg, dv_total_arr, label=r\"$\\Delta v_{total}$\", linewidth=2)\n",
"ax.axvline(angles_deg[idx_min], color=\"#BBBBBB\", linestyle=\":\", linewidth=1)\n",
"ax.set_xlabel(\"Transfer Angle (degrees)\")\n",
"ax.set_ylabel(r\"$\\Delta v$ (m/s)\")\n",
"ax.set_title(\"LEO to MEO Transfer\")\n",
"ax.legend()\n",
"plt.tight_layout()\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"id": "8",
"metadata": {},
"source": [
"## Pork-Chop Plot: Delta-v vs Transfer Angle and Time\n",
"\n",
"A **pork-chop plot** shows total delta-v as a function of both transfer angle and time of flight,\n",
"revealing the trade-off between fuel cost and transfer duration."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "9",
"metadata": {},
"outputs": [],
"source": [
"# Sweep both transfer angle and time of flight\n",
"angles = np.linspace(60, 300, 120)\n",
"tofs = np.linspace(0.4 * t_hohmann, 3.0 * t_hohmann, 100)\n",
"\n",
"DV = np.full((len(tofs), len(angles)), np.nan)\n",
"\n",
"for i, t in enumerate(tofs):\n",
" for j, theta_deg in enumerate(angles):\n",
" theta = np.radians(theta_deg)\n",
" r2 = np.array([r2_mag * np.cos(theta), r2_mag * np.sin(theta), 0.0])\n",
"\n",
" try:\n",
" solutions = sk.lambert(r1, r2, t)\n",
" v1_t, v2_t = solutions[0]\n",
"\n",
" v1_c = np.array([0.0, v_circ_1, 0.0])\n",
" v2_c = np.array([\n",
" -v_circ_2 * np.sin(theta),\n",
" v_circ_2 * np.cos(theta),\n",
" 0.0,\n",
" ])\n",
"\n",
" dv = np.linalg.norm(v1_t - v1_c) + np.linalg.norm(v2_t - v2_c)\n",
" if dv < 15000:\n",
" DV[i, j] = dv\n",
" except Exception:\n",
" pass\n",
"\n",
"fig, ax = plt.subplots(figsize=(9, 6))\n",
"levels = np.linspace(500, 8000, 20)\n",
"cs = ax.contourf(angles, tofs / 60, DV, levels=levels, cmap=\"viridis\", extend=\"max\")\n",
"ax.contour(angles, tofs / 60, DV, levels=levels, colors=\"white\", linewidths=0.3)\n",
"cbar = fig.colorbar(cs, ax=ax, label=r\"Total $\\Delta v$ (m/s)\")\n",
"ax.set_xlabel(\"Transfer Angle (degrees)\")\n",
"ax.set_ylabel(\"Time of Flight (minutes)\")\n",
"ax.set_title(\"LEO to MEO Transfer: Pork-Chop Plot\")\n",
"plt.tight_layout()\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"id": "10",
"metadata": {},
"source": [
"## Transfer Orbit Visualization\n",
"\n",
"Let's visualize the optimal transfer orbit alongside the departure and arrival orbits."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "11",
"metadata": {},
"outputs": [],
"source": [
"# Optimal transfer\n",
"theta_opt = np.radians(angles_deg[idx_min])\n",
"r2_opt = np.array([r2_mag * np.cos(theta_opt), r2_mag * np.sin(theta_opt), 0.0])\n",
"\n",
"tof_opt = t_hohmann * (theta_opt / np.pi)\n",
"solutions = sk.lambert(r1, r2_opt, tof_opt)\n",
"v1_transfer, v2_transfer = solutions[0]\n",
"\n",
"# Propagate transfer orbit via two-body integration\n",
"from scipy.integrate import solve_ivp\n",
"\n",
"def twobody(t, state):\n",
" r = state[:3]\n",
" v = state[3:]\n",
" r_norm = np.linalg.norm(r)\n",
" a = -mu / r_norm**3 * r\n",
" return np.concatenate([v, a])\n",
"\n",
"state0 = np.concatenate([r1, v1_transfer])\n",
"sol = solve_ivp(twobody, [0, tof_opt], state0, rtol=1e-10, atol=1e-12,\n",
" dense_output=True)\n",
"t_plot = np.linspace(0, tof_opt, 500)\n",
"positions = sol.sol(t_plot)[:3].T\n",
"\n",
"# Draw orbits\n",
"fig, ax = plt.subplots(figsize=(8, 8))\n",
"theta_circ = np.linspace(0, 2 * np.pi, 500)\n",
"\n",
"# Earth\n",
"earth = plt.Circle((0, 0), r_earth / 1e3, color=\"lightblue\", ec=\"steelblue\", linewidth=1)\n",
"ax.add_patch(earth)\n",
"\n",
"# Departure orbit\n",
"ax.plot(r1_mag * np.cos(theta_circ) / 1e3, r1_mag * np.sin(theta_circ) / 1e3,\n",
" \"--\", color=\"#BBBBBB\", linewidth=1, label=f\"LEO ({alt_1/1e3:.0f} km)\")\n",
"\n",
"# Arrival orbit\n",
"ax.plot(r2_mag * np.cos(theta_circ) / 1e3, r2_mag * np.sin(theta_circ) / 1e3,\n",
" \"--\", color=\"#BBBBBB\", linewidth=1, label=f\"MEO ({alt_2/1e3:.0f} km)\")\n",
"\n",
"# Transfer arc\n",
"ax.plot(positions[:, 0] / 1e3, positions[:, 1] / 1e3,\n",
" linewidth=2, color=\"#CC3311\", label=\"Transfer orbit\")\n",
"\n",
"# Endpoints\n",
"ax.plot(r1[0] / 1e3, r1[1] / 1e3, \"o\", color=\"#0077BB\", markersize=8, zorder=5)\n",
"ax.annotate(\"Departure\", xy=(r1[0] / 1e3, r1[1] / 1e3),\n",
" xytext=(15, 15), textcoords=\"offset points\", fontsize=11)\n",
"ax.plot(r2_opt[0] / 1e3, r2_opt[1] / 1e3, \"s\", color=\"#009988\", markersize=8, zorder=5)\n",
"ax.annotate(\"Arrival\", xy=(r2_opt[0] / 1e3, r2_opt[1] / 1e3),\n",
" xytext=(15, -20), textcoords=\"offset points\", fontsize=11)\n",
"\n",
"ax.set_xlabel(\"x (km)\")\n",
"ax.set_ylabel(\"y (km)\")\n",
"ax.set_title(f\"Transfer Orbit ({angles_deg[idx_min]:.0f}°, \" + r\"$\\Delta v$ = \" + f\"{dv_total_arr[idx_min]:.0f} m/s)\")\n",
"ax.set_aspect(\"equal\")\n",
"ax.legend(loc=\"lower left\")\n",
"plt.tight_layout()\n",
"plt.show()"
]
},
{
"cell_type": "markdown",
"id": "12",
"metadata": {},
"source": [
"## Interplanetary Example: Earth to Mars\n",
"\n",
"Lambert targeting also works for interplanetary transfers by using the Sun's gravitational parameter and heliocentric positions."
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "13",
"metadata": {},
"outputs": [],
"source": [
"mu_sun = sk.consts.mu_sun\n",
"au = sk.consts.au\n",
"\n",
"# Simplified circular coplanar orbits\n",
"r_earth_orbit = 1.0 * au\n",
"r_mars_orbit = 1.524 * au\n",
"\n",
"# Earth at departure\n",
"r1_helio = np.array([r_earth_orbit, 0, 0])\n",
"\n",
"# Mars at arrival (~135 degree transfer)\n",
"transfer_angle = np.radians(135)\n",
"r2_helio = np.array([\n",
" r_mars_orbit * np.cos(transfer_angle),\n",
" r_mars_orbit * np.sin(transfer_angle),\n",
" 0,\n",
"])\n",
"\n",
"# Transfer time: ~8 months\n",
"tof_days = 245\n",
"tof_interp = tof_days * 86400\n",
"\n",
"solutions = sk.lambert(r1_helio, r2_helio, tof_interp, mu=mu_sun)\n",
"v1_t, v2_t = solutions[0]\n",
"\n",
"# Circular velocities\n",
"v_earth = np.sqrt(mu_sun / r_earth_orbit)\n",
"v_mars = np.sqrt(mu_sun / r_mars_orbit)\n",
"\n",
"v1_earth = np.array([0, v_earth, 0])\n",
"v2_mars = np.array([\n",
" -v_mars * np.sin(transfer_angle),\n",
" v_mars * np.cos(transfer_angle),\n",
" 0,\n",
"])\n",
"\n",
"dv_depart = np.linalg.norm(v1_t - v1_earth)\n",
"dv_arrive = np.linalg.norm(v2_t - v2_mars)\n",
"\n",
"print(f\"Earth-Mars transfer ({tof_days} days):\")\n",
"print(f\" Departure delta-v: {dv_depart/1e3:.2f} km/s\")\n",
"print(f\" Arrival delta-v: {dv_arrive/1e3:.2f} km/s\")\n",
"print(f\" Total delta-v: {(dv_depart + dv_arrive)/1e3:.2f} km/s\")"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"name": "python",
"version": "3.12.0"
}
},
"nbformat": 4,
"nbformat_minor": 5
}